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Abstract. We present a Lagrange-Galerkin scheme free from numerical quad¬ 
rature for convection-diffusion problems. Since the scheme can be implemented 
exactly as it is, theoretical stability result is assured. While conventional 
Lagrange-Galerkin schemes may encounter the instability caused by numer¬ 
ical quadrature error, the present scheme is genuinely stable. For the P^- 
element we prove error estimates of 0(At-\-h^ in -norm and of 

0(At -{- -{- h^) in £°°(iL^)-norm. Numerical results reflect these estimates. 


1. Introduction 

The Lagrange-Galerkin method, which is also called characteristics finite ele¬ 
ment method or Galerkin-characteristics method, is a powerful numerical method 
for flow problems such as the convection-diffusion equations and the Navier-Stokes 
equations. In this method the material derivative is discretized along the charac¬ 
teristic curve, which originates the robustness for convection-dominated problems. 
Although, as a result of the discretization along the characteristic curve, a compos¬ 
ite function term at the previous time appears, it is converted to the right-hand 
side in the system of the linear equations. Thus, the coefficient matrix in the left- 
hand side is symmetric, which allows us to use efficient linear solvers for symmetric 
matrices such as the conjugate gradient method and the minimal residual method 

mus]. 
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Stability and error analysis of LG schemes has been done in [D n n El m iini 
nn m m m HU; see also the bibliography therein. Pironneau m analyzed 
convection-diffusion problems and the Navier-Stokes equations to obtain subopti- 
mal convergence results. Optimal convergence results were obtained by Douglas- 
Russell [6] for convection-diffusion problems and by Siili m for the Navier-Stokes 
equations. Optimal convergence results of second order in time were obtained by 
Boukir et al. [4] for the Navier-Stokes equations in multi-step method and by 
Rui-Tabata m for convection-diffusion problems in single-step method. All these 
theoretical results are derived under the condition that the integration of the com¬ 
posite function term is computed exactly. Since, in real problems, it is difficult 
to get the exact integration value, numerical quadrature is usually employed. It 
is, however, reported that instability may occur caused by numerical quadrature 
error in [S1II71IIH]. That is, the theoretical stability results may collapse by the 
introduction of numerical quadrature. 

Several methods have been studied to avoid the instability. The map of a particle 
from a time to the previous time along the trajectory, which is nothing but to 
solve a system of ordinary differential equations (DDEs), is simplified in P fT^. 
Morton-Priestley-Suli [9] solved the ODEs only at the centroids of the elements, 
and Priestley m did only at the vertices of the elements. The map of the other 
points is approximated by linear interpolation of those values. It becomes possible 
to perform the exact integration of the composite function term with the simplified 
map. Bermejo-Saavedra [3] used the same simplified map as [T3] to employ a 
numerical quadrature of high accuracy to the composite function term. Tanaka- 
Suzuki-Tabata m approximated the map by a locally linearized velocity and the 
backward Euler approximation for the solution of the ODEs in Pi-element. The 
approximate map makes possible the exact integration of the composite function 
term with the map. Pironneau-Tabata [12] used mass lumping in Pi-element to 
develop a scheme free from quadrature for convection-diffusion problems. 

In this paper we prove the stability and convergence for the scheme with the same 
approximate map as m in Pfc-element for convection-diffusion problems. Since we 
neither solve the ODEs nor use numerical quadrature, our scheme can be precisely 
implemented to realize the theoretical results. It is, therefore, a genuinely stable 
Lagrange-Galerkin scheme. Our convergence results are of 0{At + h'^ + in 

^oo(2^2)-noi-m and of 0{At + h'^ + h^) in £°“(iJ^)-norm. They are best possible in 
both norms for Pi-element and in (iJ^ )-norm for P 2 -element 

The contents of this paper are as follows. In the next section we describe the 
convection-diffusion problem and some preparation. In section |3l after recalling the 
conventional Lagrange-Galerkin scheme, we present our genuinely stable Lagrange- 
Galerkin scheme. In section|4]we show stability and convergence results, which are 
proved in section [5l In section [6] we show some numerical results, which reflect the 
theoretical convergence order. In section [3 we give conclusions. 

2. Preliminaries 

We state the problem and prepare notation used throughout this paper. 

Let D be a polygonal or polyhedral domain of (d = 2, 3) and T > 0 be a time. 
We use the Sobolev spaces LP{Vt) with the norm IHlop, IT®’P(n) and ITg’^(D) with 
the norm IHlgp and the semi-norm 1-1^^ for 1 < p < oo and a positive integer s. 
We will write H^{n) = and drop the subscript p = 2 in the corresponding 
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norms. The L^-norm H-Hg is simply denoted by H-H. The dual space of is 

denoted by For the vector-valued function w G define the 

semi-norm |w|i^oo by 



0,oo 


The parenthesis (•,•) shows the L^-inner product (/, g) = /n/ff ^or a Sobolev 
space X(n) we use abbreviations H™‘{X) = iJ™(0, T; X(n)) and C{X) = (^([O, T]]X{n)). 
We define a function space Z'^{ti,t 2 ) by 


Z^ihM) = {/e = 0,.. 

2 ■ 

''j=0 




and denote Z^{0,T) by Z^. 

We consider the convection-diffusion problem: find : fl x (0, T) —>■ R such that 


^+u-V(p-vA(p = /, 

(x, t) G 12 X (0, T), 

(la) 

(P = 0, 

(x,t) G ai2 X (0,T), 

(lb) 

o 

II 

X G 12, t = 0, 

(Ic) 

where is the boundary of 12 and > 0 is 

a diffusion constant which is 

less than 


or equal to a given uq. Functions u : fl x (0,T) —>• / G C{L'^) and G (7(^2) 

are given. 

Remark 1. As usual, in place of (HEl, we can deal with the inhomogeneous boundary 
condition cj) = g hy replacing the unknown function (j)hy (f> = (f> — g ii the function g 
dehned on x (0,T) can be extended to a function g in x (0,T) appropriately. 


Let At > 0 be a time increment, Nt = [T/AtJ, t" = nAt and for 

a function ip defined in ft x (0, T). For a set of functions ip = {V'"'}^=o’ norms 

and are defined by 

||V’llf=o(i 2 ) = max{||i/’”|| ;n = 0 ,..., At}, 

C n2 

At E ll^”ll' 

n—ni 


and denote ||t^||^ 2 (i_Ar,^.i 2 ) by \\ip\\e 2 (^L^y 

Let u be smooth. The characteristic curve X(t;x,s) is defined by the solution 
of the system of the ordinary differential equations, 

^^(t; cc, s) = ulAlt; X, s), t), t < s, (2a) 

at 

A(s; X, s) = X. (2b) 


Then, we can write the material derivative term ^ + u - X(p as 
(^+«-V0) (X(t),t) = |<^(A(t),t). 
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For ui : ri —>■ we define the mapping Xi('u;) : —>■ by 

(Xi{w))(x) = X — w(x)At. (3) 

Remark 2. The image of x by Xi(u{-,t)) is nothing but the backward Euler ap¬ 
proximation oi X{t — At; X, t). 

The symbol o stands for the composition of functions, e.g., {go f){x) = g{f{x)). 
Let 7/t = {K} be a triangulation of and h = maxKeTh diam(iL) be the maxi¬ 
mum element size. Throughout this paper we consider a regular family of triangu¬ 
lations {Th}hio- Let fc be a fixed positive integer and 14 C HqIH) be the Pfc-finite 
element space, 

Vh = {Vh G c{n) n H^{n); v^k g Pk{K), vk g r4, 

where Pk{K) is the set of polynomials on K whose degrees are less than or equal 
to k. Let (j)h G 14 be the Poisson projection of </) G defined by 

^^h€Vh. (4) 

We use c to represent a generic positive constant independent of h, At, v, f and 
(j) which may take different values at different places. The notation c{A) means that 
c depends on a positive parameter A and that c increases monotonically when A in¬ 
creases. The constants Co, ci and C 2 stand for cq = c(||u||(;i(-j;^oo)), ci = c(||u||p(-^yi,cx,)) 
and C 2 = c{\\u\\(yf^^ 2 ,oo'j)- We also use fixed positive constants a* and defined in 
Lemma [T] in the next section and in Lemma [S] in Section [SJ respectively. 

3. A GENUINELY STABLE LAGRANGE-GALERKIN SCHEME 

The conventional Lagrange-Galerkin scheme, which we call Scheme LG, is de¬ 
scribed as follows. 

Scheme LG. Let = 0° . Find C 14 such that for n = 1,..., Nt 

o X" \ 

- = ViA, GP^, (5) 

where Af = Ai(u"). 

For this scheme error estimates 

Uh - < c{h’^ + At), c(l/u)(/i'=+i -h At), 

Uh - (t>\\i^(m) < c(l/u)(/i'= -f At) 

are proved in [^, where the composite function term {4>^~^ o A", 'iph) is assumed to 
be exactly integrated. 

Although the function ^ polynomial on each element K, the composite 

function o A" is not a polynomial on K in general since the image A"(A) 
of an element K may spread over plural elements. Hence, it is hard to calculate 
the composite function term o A”, 4ft) exactly. In practice, the following 

numerical quadrature has been used. Let (7 : A —>■ R be a continuous function. A 
numerical quadrature Ih[g; K] of gdx is defined by 

Ih[g;K] =meas(A)y^Wi gja^), (7) 
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where Nq is the number of quadrature points and (wi,ai) G K. x if is a pair of 
weight and point for i = 1,..., Nq. We call the practical scheme using numerical 
quadrature Scheme LG'. 

Scheme LG'. Let = 0°. Find C 14 such that for n = 1,..., Nt 

Ken ( 8 ) 

= {r,^Ph), V4/.eL/., 


where = Xi(u"). 


It is reported that the results ([5]) do not hold for Scheme LG' [SJ HZl UHl US] ■ 

We denote by the Lagrange interpolation operator to the Pi-finite element 
space. The following lemma is well-known [^. 

Lemma 1. (i) There exists apositive eonstantcu such thatforw G 

- w||o,oo < cn/i^ kla.oo ■ 

(ii) There exists a positive constant a* > 1 such that for w G 
|nj,^^w|i,oo < a* . 

We now present our genuinely stable scheme GSLG, which is free from quad¬ 
rature and exactly computable. We define a locally linearized velocity Uh and a 
mapping by 

Scheme GSLG. Let (jf^ = 4>°. Find {4>h}n=i for n = 1,..., Nt 

( An _ X"—1 o X" \ 

= V4/.G14. (9) 

We show that the integration (4t~^ ° can be calculated exactly. At 



Figure 1. Elements Kq, Ki and a polygon Ei 

first we prepare two lemmas. The next lemma on the mapping (jS)) is proved in |14) . 

Lemma 2 ([141 Proposition 1]). Suppose 

w G and At oo 1- (^9) 

Let F = Ai(rt;) be the mapping defined in (0). Then, F : LI ^ Ll is bijective. 

Lemma 3. Let Kq, Ki G 7/t and F : Kq —>■ R'* be linear and one-to-one. Let 
Fi = Kq n F~^{Ki) and meas(Ei) > 0. Then, the following hold. 
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(i) El is a polygon (d = 2) or a polyhedron (d = 3). 

(ii) (j)h o G Pk{Ei), V4>h G Pk{Ki). 

Proof, (i) Since both Kq and F~^{Ki) are triangles {d = 2) or tetrahedra (d = 3), 
the intersection is a polygon or a polyhedron. See Fig. [T] 

(ii) F G Pi{KqY implies that F G Pi{EiY and it holds that F{Ei) C Ki. 
Hence, (fih o F^^i is well defined and (j)h o G Pk{Ei). □ 

Proposition 1. Let (f>h, V'/i G ^h, w G Wg^’°“(r2) and Xih = where Xi 

is the operator defined in (0). Suppose a*At < 1- Then, o Xih)'iphdx 

is exactly computable. 


Proof. It is sufficient to show that ° Xih)'tphdx can be computable exactly 

for any Kq G Tj). The mapping Xi^ : H —>■ H is bijective since we can apply Lemma 
[5] thanks to 

At|n|^^^w|i,oo < a^,At\w\i^oo < 1- (11) 

Let A(Fo) = {/; Ko n Xfi^\Ki) ^ 0} and F, = KodXfi^^Ki) for I G A{Ko). Noting 
that 

U El = Ron U Xfi^\Ki) = Ko 

ieA(ifo) /GA(ifo) 

and that meas(F; fl Em) = 0 ioi I m, l,m G A(Ko), we can divide the integration 
on Kq into the sum of those on Ei for I G A{Kq), 


’Kq 


{4>h o Xih)i>hdx = 


/ {fihO Xik)iphdx. 

leA(Ko) 


Since Lemma [3] with F = Xik implies that both fik ° Xih and ifh are polynomials 
on El, we can execute the exact integration. □ 


Remark 3. In the case of d = 2, Priestley [T3] approximated A(F a;,t") by 

Xih{x) = BiXi{x) + B2X2 {x) + B^Xzix), X G Kq 

on each Kq G Th, where Bi = A(f"“^;Hi,t”), are vertices of Kq and 

{Ai}^^i are the barycentric coordinates of Kq with respect to Since Xih{x) 

is linear in Kq, the decomposition 

/ {(j)h o Xih)tphdx = / {fiho Xihfi’hdx, 

•^^0 leA(Ko) 

A(Fo) = {/; Kq n Xfi^\Ki) ^ 0} , F, = Fq n Xfi^\Ki) 

makes the exact integration possible. However, Bi = A(t”“^; f") are the solu¬ 

tions of a system of ordinary differential equations and they cannot be solved ex¬ 
actly in general. In practice, some numerical method, e.g., Runge-Kutta method, 
is required, which introduces another error. 


4. Main results 


We show the main results, the stability and convergence of Scheme GSLG. 
Hypothesis 1. (i) u G C{{Wl^^)^), (ii) u G C{{Wl^^ n 1^2’°°)^). 

Hypothesis 2. 0 G n 
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Hypothesis 3. The time increment At satisfies 0 < At < Ato, where 

Ato = ^-, ( 12 ) 

oi* mc(w^’°°) 

and a* and d* are the constants stated in Lemma [T] (Section [3|) and Lemma [5] 
(Section [S]), respectively. 

Hypothesis 4. There exists a positive constant cp such that, for ip G fl 

Wi’h - V'llo < llV'll /c+1 7 (13) 

where iph is the Poisson projection defined in (|^. 

Remark 4. (i) It is well-known that the iL^-estimate 

- V'lli < ||^/>||^_^^ (14) 

holds without any specific condition. On the other hand, Hypothesis H] holds, for 
example, if H is convex, by Aubin-Nitsche lemma 
(ii) Hypothesis [5] implies (p G and <p^ G 

Theorem 1. Suppose Hypotheses{^(i) and\^ Then, there exists a positive constant 
c* independent of h, At, v, (p and f such that 


ll'(’?*ll^“(L2) + Vv\\y (ph\\l2(^p2) < c* ^||()>°|| + ||/||^2(j;,2)^ . 

Theorem 2. Suppose Hypotheses[J\-(ii) ,1^ andl^ 

(i) There exists a positive constant c* independent of h,At,v and (p such that 

\\<P — 4>h\\£ac(^p2^ + \fv ||V((^ — Ph)\\i2(^p2'^ 

dp 


at 


<c*| At \\p\\z 2 +h'^(^ 


+ h? ||V(/)||^ 2 (o, 7 VT-l;i=*)|- 


P 2 (pfe + 1 ) 


+ 11011 + II0II 

£2(o,7VT;ff''+i)) (15j 


(a) There exists a positive constant c** independent of h,At,p (but dependent on 
\/v) such that 


( 

dp 

V 

at 


110 Ph\\)ioo(pi\ + At II 0 II 22 + h 

I \ UL p2(^fc + l) 

+ II0II \ G-hf ||V0||^2(o_7Vt-1;L2) f ■ 


+ II0II 


^“(P'=+ 1 ) 


(16) 


(Hi) Moreover, suppose Hypothesis^ Then, there exists a positive constant c**, 
independent of h, At, v and p such that 


110 0ft||foo(j;,2) + At 11011^2-t ^ 


dp 


at 


L 2 (i?'= + l) 


+ ll 0 L“(p'=+i) 


+ ^-'/^ll 0 ll.. 


(17) 
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Remark 5. From Theorem!^ we have 


\\(l> - < c{At + h'^ + h^), c (^At + h"^ + 

||(/) - 0?i||^=c.(^i) < c {At + h'^ + h^). 

In the case of Pfc-element, k = 1,2, the estimate (USD shows the optimal L^- 
convergence rate 0{At + h^) independent of v. The dependency on v in m 
and (ED is also inevitable in Scheme LG. 


5. Proofs of main theorems 


We recall some results used in proving main theorems. For their proofs we only 
show outlines or refer to the bibliography. 

Lemma 4 ([TH Lemma 1]). Suppose w € and 

(18) 

Let F = Xi {w) be the mapping defined in (0). Then, there exists a positive constant 
c(|r(;|ioo) such that for tp £ (LI) 

Wi’oFW < (l+cAt) ||i/;||. 


The proof is given in El- 

Lemma 5. There exists a constant S^, £ (0,1) such that, for w £ lFQ^’°°(n)'^ and 
At satisfying At |w| ^ < <5*, 


1 ^ dXi{w) 

2 ~ dx 



where \dXi{w)/dx\ is the Jacobian of the mapping Xi{w) defined in (0)- 


Lemma [5] is easily proved by the fact, 


/ dXi{w) 
V dx 


- Q ■ 


Lemma 6. Let Wi £ Wg^’°°(fl)‘^ and F^ = Xi{wi) be the mapping defined in 13) 
for i = 1,2. Under the condition At oo — ~ 1’ have for tp £ 


Ilf^oFi -f^o F 2 II < V^At llrci - W2||o,oo II^V'll ■ 

Lemma ini is a direct consequence of [TJ Lemma 4.5] and Lemma [SJ 

Lemma 7. Let w £ lFQ^’°“(n)‘^ and F = Ai(ui) be the mapping defined in 
Under the condition < (5*, there exists a positive constant c(||r(;||^^) 

such that for ip £ L‘^{Lt) 


\\iP-iPoF\\h~^qj < cAt\\ip\\. 


Lemma [7] is obtained from (BJ Lemma 1] and Lemma [S] 






Lemma 8 (discrete Gronwall inequality). Let oq and oi be non-negative numbers, 
At G (0, be a real number, and {a:”}„>g , {y”}„>i and be non-negative 

sequences. Suppose 

™n ™n—1 


At 


+ y” < aox^ + + 6”, Vn > 1. 


Then, it holds that 


n / n 

< exp{(2ao + ai)nAt} ( + At ^ 6* 

i=i \ i=i 


Vn > 1. 


Lemma |5] is shown by using the inequalities 
1 


< 1 + 2aoAt < exp(2aoAt). 


1 — ttoAt 

Outline of the proof of Theorem [7] We substitute (j)^ into ifh in We can apply 
LemmalHwith w = and by virtue of At < 1- The rest of 

the proof is similar to [HI Theorem 1]. We, therefore, omit it. 

Proof of Theorem\^ We hrst show the estimate (fT51) . Let 

ah — f'h V — 4^ 4*h: 

where (ph is the Poisson projection defined in (|T]). From © and m we have 


n_ n-l y 
Cu e. u yv- 


Ih 


At 


for iph & Vh, where 


dp” 


R” = -^ + u”- yp” - 
at 


Ph]+i^{^el,y^fh) = Y.{R7.P’h) 

1=1 

-P”-^ o X” 


( 20 ) 


At 


on — 

ri 2 — 


R” = 


Ip- - T] 


oX”^-P”-^oXf 

At 

_ r^ri-1 _ r^n-1 o xn^ 


( 21 ) 


on — 
-tCa = 


At ’ ^ “ At 

Substituting into 'pn, applying Lemma |4] with F = and p = epp^, and 
evaluating the first term of the left-hand side as 

n—1 

:) /A . \ \ I .. T ... 

— n X” 1 

\Ph ° ^lh\\ ) 


p“' _ p"' Q 

^ ^Ih 


At 


>^(l|e)(r-(l + ciAt)^||e;: 


1 


-)-^(2 + ciAt)||e;; 


we have 


2At 


<ci e" 


eh - Ch 


■i^\\ye”h 




(22) 


2=1 


^2 = 1 


where {ei}i=i are positive constants satisfying Ato < Eq = Hi=i 










We evaluate i = 1, • • ■ ,4. Setting 


y{x,s) = X + {s — l)At u'^{x), t{s)=t^ ^ + sAt, 


we have 


which implies 


/>" - o 1 


At 




s=0’ 


d(j}” 


i?? = ^ + u” • - 

at 


u”(-) • V(?!) + ^ i- (y(-, s), t(s))(is 


(9t 


= At ds 





= Aty I (2/(-,si),i(si))dsi. 


Hence, we have 


||i?5^|| < AtJ^' si I (^u"(-) . V + (y(.,si),t(si)) 


dsi 


< Co A/At ||(^|| 




(23) 


where we have used the transformation of independent variables from x to y and 
Si to t, and the estimate \dx/dy\ < 2 by virtue of Lemma [S] 

From At |u|(^^^i,oo) , At\uh\(j(yY^’°°) — ^*’ Lemmas |6] and [1] it holds that 

||i?^|| < y2||V</)"-i||||n«u" -u"||o.oo < C 2 h^ . (24) 


i ?3 is evaluated as 


Pall 



dij 


{■,t{s))ds 


< 

“ ^At 


d4> 

'm 


L2(tn-ytn.^ffc + l) 


where we have used m- 

From At |rt/i|c(vi^i,cx,) < and Lemma[6]it holds that 


p4”ll < y2||Vr?'^-i||nWu"||o,oo < coh'= 
Combining ([2^ - (l26ll . we have 


^ {\Kf - ||eri^) + ^ II Vejjf < £0 \\elf + c, \\el 


+ 


'2k 


At 


d(j) 

'dt 


2 




2 

fe+1 


(25) 


(26) 
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From Lemma [5] we obtain for n = 1,..., Nt 

2 


< C 2 




Nt 

+ 2 i/At^ 
i=i 

, 2 AfT -1 

at 

which implies (1151) by virtue of = 0 and the triangle inequalities, 


A^t —1 

+ h^At Y, \\y<jyi\f 

j=o 


||<^ — 4>h\\eoc(^i^2) < \\eh\\i^(^L^') + lhll^=x=(L2) 

< ||eft,||^oo(J^2) + Cp/l^ ||</'||^oo(pfc + l) , 

||V(^ — ^h)||^2(p2) < II Ve?i||^2(p2) + ||V?7||^2(p2) 

< IIVeh||^2(p2) + cph^ ll^ll^2(p-fc+i) • 

We show the estimate (fT 6 l) . The equation (l 20 l) can be rewritten as 

1 5 


(27) 


— {el - el-\i,H) + KVe)), V^-^) = ^h), 


where 




From Lemma [5] it holds that 

||i?^|| < v^||Veri||||n(V||o.oo < co||Veri||. 

Substituting D^tel = '^i^l ~ ^h~^) V'/u and using (l^ - (05l) for Ri,..R 4 , 
we have 


- - Ve 


||i?A*e)(f+ ^Q||Ve 

<C2 |At ||0|||2(,„-y,„) + + 


V(e)(- 


n n—l\ 


At 


d(j) 


2At 

2 


dt 




Ve' 


n — 1 


+'-“lk"-‘llL.} + 7( 

From Lemma [5] we have for n = 1, -. -, A^t 

-T Y II^Ate)(f + ^ ||Ve)(f < C2 exp 




i=i 


(f 


Vel\\^ + AfM% 


+ h 


2k 


d(f) 


dt 


Wt —1 iVT-- —1 -V 

+ k“At y; 11^11;^, +ft-At y; iiv^in, 

i=o j=o ^ 


— 1 


p2(pfe + l) 

which implies (1161) by virtue of e® = 0 , the triangle inequality, 

||V((^ — (^?i)||^oo(p2) < IIVe/i||^cxj(p2) + IIVr7||^oo(p2) 

< IIVe/i||^oo(p2) + cph^ ll'i^ll^=x=(p-fc+i) 

and the Poincare inequality, 

|lu||i<c||Vn||, yvGH^in). (28) 
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Now we show the estimate CZD- We return to the error equation (I^Dll . Using 
m in place of (fni) in the estimate of K^, we can evaluate (05|) as 

d(j) 


|i?3 II < 


\fKi 


dt 




From Lemma [7] we have 


\\R2\\H-Hn) < Cl < cih>^+^ \\r-^\k+i ■ 

Hence, it holds that 

{R2,el) < ||i?4”ll^-i(n) WelWi < \\r-^\\l+i + ^ llVe^jf, 

where we have used the Poincare inequality (1281) . Using this inequality instead of 
^ ||i ?4 ||^+e 4 ||e}J||^ in (|2^ and replacing the last term of (|27l) by ll?^ll^oo(fffc+i)) 

we obtain dm). 

□ 


6. Numerical results 


We show numerical results in d = 2. We compare the conventional scheme 
(Scheme LG') with the present one (Scheme GSLG). We use FreeFem++ [8] for 
the triangulation of the domain. Both Pi- and P 2 -elements are used. For Scheme 
LG' we use the seven points quadrature formula of degree five [7]. A relative error 
Ex is defined by 

||n("Vll^»(x) 

(k) 

where n)^ is the Lagrange interpolation operator to the Pfe-finite element space 
and X = Lp'{Vt) or Pq(H). 


Example 1 (The rotating Gaussian hill [M]). In ([T|), H is an unit disk, and we set 
T = 2tt,v = 10 -^ 

u{x,t) = {-X2,xi), f = 0, (l)° = (j)e{-,0), 


where 


J. , P cr f (xi(t) - Xi^cf-h (X2(t) - X2,cy 

9 e{x, t) = —exp -——- 

a + ivt ( cr + 

{xi , X 2 ) (t) = {xi cos t X 2 sin t,—xi sin t X 2 cos t ), 

(iCi.c, X2,c) = (0.25, 0), a = 0.01. 


In this problem the identity hI^^^u = u holds. This problem does not satisfy 
our setting because H is not a polygon and u ^ 0 on dfl. The function (pe in Fig. 
[7] (left) satisfies (TTal) and dS but does not satisfy the boundary condition (HEl. 
However, we may apply the schemes and treat (pf, as the solution since the value of 
(j>e on is almost equal to zero, less than 10“^^, and we may neglect the effect of 
the boundary value and the term o Xi)iphdx and o X^f^)ip/idx on 

the element K touching the boundary. 

Let N be the division number of the circle. We set h = 27r/N, N = 32,64,128 and 
256. Figured] (right) shows the triangulation of H for = 64. The time increment 
At is set to be cih and C 2 h'^ for Pi-element (ci = ^ '= 0.255, C 2 = =. 1.30), 
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Figure 2. The function (/)e(-,0) (left) and the triangulation of 
for TV = 64 (right) in Example [T] 

csh^ and for P 2 -element (ca = = 2.59, C 4 = = 13.21) so that we can 

observe the convergence behavior of 0 {h^) for and 0 {h^) and 0 (/i^+^) for 

Ei ,2 when P^-element is employed. 

In the following figures we use the symbols shown in Table 1. Figure [3] shows 
the log-log graphs of and E^i versus h. The left graph shows the results of 
Pi-element and Tables [2] shows the values of them. The convergence order of Ei ,2 
with At = 0{h) is less than 1 in Scheme LG' (C) and more than 1 in Scheme 
GSLG ( ). The orders of Ei ^2 with At = 0{h?') are almost 2 for small h in both 
schemes (E, ). The convergence of E^i is not observed in Scheme LG' (■) while 

the order is almost 1 in Scheme GSLG (•). The right graph of Fig. [3] shows the 
results of P 2 -element and Tables [3] shows the values of them. The errors Ei ^2 with 
At = 0{h?) are too large at A = 128 and 256 to be plotted in the graph in Scheme 
LG' (C) while the convergence order is almost 2 in Scheme GSLG ( ). The error 
Ei ^2 with At = 0{h^) is large at N = 128, but it becomes small again at N = 256 
in Scheme LG' ( 0 ). We will discuss the reason why such a behavior occurs in a 
forthcoming paper. The order is greater than 2.5 in Scheme GSLG ( ). The errors 
E^i are too large at A = 128 and 256 to be plotted in the graph in Scheme LG' (■) 
while we can observe the convergence of E^ji but the order is less than 2 in Scheme 
GSLG (•). The errors of Scheme GSLG are smaller than those of Scheme LG' in 
both cases of Pi- and P 2 -element. 

Table 1 . Symbols used in Figs. |3l |6] and [7] and Tables [2HZl 


A 




At 

0 {h>^) 

0 (/r'=+i) 

0 {h^) 

Scheme LG' 

□ 

1 ] 

■ 

Scheme GSLG 



• 


Figure m shows the solutions (j)^ for h = 27r/64 = 0.0982, At = 0.0065, nAt '= 2?!. 
In the case of Pi-element, the solution of Scheme LG' is oscillatory while that of 
Scheme GSLG is much better though a small ruggedness is observed. In the case 
of P 2 -element, the solution of Scheme LG' is quite oscillatory while that of Scheme 
GSLG is stable. 
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h h 

Figure 3. Graphs of and E^ji versus h in Example [1] by 
Pfc-element. k = 1 (left) and k = 2 (right) 




Figure 4. Solutions (j)^ (riAt '= 2 tt) in Example [I] by Scheme 
LG'(top left) and Scheme GSLG(top right) for Pi-element, and by 
Scheme LG'(bottom left) and Scheme GSLG(bottom right) for P 2 - 
element 
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Table 2. The values of errors and orders of the graph in Fig. |3] 
by Pi-element 


N 

□ 

order 

E 

order 

■ 

order 

32 

7.58E-01 


7.58E-01 


9.52E-01 


64 

5.65E-01 

0.42 

5.53E-01 

0.45 

1.04E-P00 

-0.13 

128 

3.93E-01 

0.52 

1.87E-01 

1.56 

9.72E-01 

0.10 

256 

2.04E-01 

0.95 

4.15E-02 

2.17 

5.54E-01 

0.81 

N 

C; 

order 


order 

• 

order 

32 

7.26E-01 


7.26E-01 


9.01E-01 


64 

4.28E-01 

0.76 

4.18E-01 

0.80 

7.34E-01 

0.30 

128 

1.36E-01 

1.65 

1.04E-01 

2.01 

3.07E-01 

1.26 

256 

5.62E-02 

1.27 

2.43E-02 

2.10 

1.45E-01 

1.08 


Table 3. The values of errors and orders of the graph in Fig. [3] 
by P 2 -element 


N 

□ 

order 


order 

■ 

order 

32 

6.86E-01 


6.86E-01 


9.22E-01 


64 

4.06E-01 

0.76 

3.97E-01 

0.79 

1.31E-P00 

-0.51 

128 

1.67E-^02 

-8.68 

8.30E-01 

-1.06 

1.72E-P03 

-10.36 

256 

1.42E-H27 

-82.81 

3.05E-03 

8.09 

3.10E-P28 

-83.90 

N 

C 

order 

© 

order 

• 

order 

32 

6.03E-01 


6.03E-01 


7.38E-01 


64 

2.09E-01 

1.53 

1.17E-01 

2.37 

3.33E-01 

1.15 

128 

5.48E-02 

1.93 

1.59E-02 

2.88 

9.86E-02 

1.76 

256 

1.38E-02 

1.99 

2.66E-03 

2.58 

3.97E-02 

1.31 


Example 2. In ([T]), ft is the square (0,1) x (0,1), and we set T = 1, z/ = 10 ^ and 
10 -^ 


u{x,t) = (sinTrai sin7ra;2,sin7ra;i sin7ra:2), / 


d(j)e _ , 

— -f U ■ 

at 




where 

(j)e{x, t) = cos(27rt) sin^(7ra:i) sin(27rx2). 

In this problem, 7 ^ u. Let N be the division number of each side of fi. 

We set h = 1/N,N = 8,16,32 and 64. Figure [5] shows the triangulation of fl 
for N = 16. The time increment At is set to be cih and 02 / 1 ^ for Pi-element 
(ci = 0.125, C 2 = 1), czh^ and for P 2 -element (03 = 1 ,C 4 = 5.12) so that we 
can observe the convergence behavior of 0 {h^) for Ej^i, and 0 {h^) and 
for E ]^2 when Pfe-element is employed. 

Figure [5] shows the log-log graphs of Ei ^2 and Efji versus h with v = 10“^. The 
left graph shows the results of Pi-element and Table U] shows the values of them. 
The convergence orders of E ]^2 with At = 0{h) are almost 1 in both schemes (□, 
. ). The orders of Ei ^2 with At = 0{h'^) are almost 2 in both schemes (0, S). The 
orders of E^i are almost 1 in both schemes (■, •). The right graph of Fig. | 6 ] shows 
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Figure 5. The triangulation of Q ior N = 16 in Example [2] 




h h 

Figure 6 . Graphs of £'^2 and Ej^i versus h in Example H] with 
1 / = 10“^ by Pfc-element. k = 1 (left) and k = 2 (right) 




h h 

Figure 7. Graphs of £^2 and Ej^i versus h in Example [5] with 
V = 10“® by £fe-element. k = 1 (left) and k = 2 (right) 
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Table 4. The values of errors and orders of the graph in Fig. |6] 
by Pi-element 


N 

□ 

order 

0 

order 

■ 

order 

8 

8.14E-02 


8.14E-02 


l.lOE-01 


16 

3.64E-02 

1.16 

1.90E-02 

2.10 

4.36E-02 

1.34 

32 

1.70E-02 

1.10 

4.58E-03 

2.05 

1.87E-02 

1.22 

64 

8.53E-03 

0.99 

1.19E-03 

1.94 

9.18E-03 

1.03 

N 


order 


order 

• 

order 

8 

8.97E-02 


8.97E-02 


1.09E-01 


16 

3.68E-02 

1.29 

2.13E-02 

2.07 

4.23E-02 

1.37 

32 

1.78E-02 

1.05 

4.83E-03 

2.14 

1.92E-02 

1.14 

64 

8.90E-03 

1.00 

1.29E-03 

1.90 

9.43E-03 

1.03 

Table 5. The values of errors and orders of the graph in 

Fig. m 

by P 2 -element 






N 

□ 

order 

0 

order 

■ 

order 

8 

7.01E-02 


4.50E-02 


7.37E-02 


16 

1.77E-02 

1.99 

5.72E-03 

2.98 

1.85E-02 

1.99 

32 

4.44E-03 

2.00 

7.11E-04 

3.01 

4.63E-03 

2.00 

64 

l.llE-03 

2.00 

8.86E-05 

3.00 

1.15E-03 

2.01 

N 


order 


order 

• 

order 

8 

6.31E-02 


3.87E-02 


6.60E-02 


16 

1.58E-02 

2.00 

6.89E-03 

2.49 

1.64E-02 

2.01 

32 

3.98E-03 

1.99 

1.42E-03 

2.28 

4.12E-03 

1.99 

64 

9.90E-04 

2.01 

3.37E-04 

2.08 

1.02E-03 

2.01 

Table 6 . The values of errors and orders of the graph in 

Fig. m 

by Pi-element 






N 

□ 

order 

0 

order 

■ 

order 

8 

9.53E-02 


9.53E-02 


2.90E-01 


16 

3.94E-02 

1.27 

4.17E-02 

1.19 

2.00E-01 

0.54 

32 

1.82E-02 

1.11 

1.33E-02 

1.65 

1.09E-01 

0.88 

64 

9.07E-03 

1.00 

3.44E-03 

1.95 

5.12E-02 

1.09 

N 


order 


order 

• 

order 

8 

9.70E-02 


9.70E-02 


2.23E-01 


16 

3.93E-02 

1.30 

2.86E-02 

1.76 

1.23E-01 

0.86 

32 

1.89E-02 

1.06 

7.42E-03 

1.95 

6.02E-02 

1.03 

64 

9.45E-03 

1.00 

1.91E-03 

1.96 

3.08E-02 

0.97 


the results of P 2 -element and Table [5] shows the values of them. The convergence 
orders of Pl 2 with At = 0{h?) are almost 2 in both schemes (□, ). The order 

of Ei ^2 with At = 0{h^) is almost 3 in Scheme LG' (0) and almost 2 in Scheme 
GSLG ( ). The orders of Ejji are almost 2 in both schemes (■, •). These results 
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Table 7. The values of errors and orders of the graph in Fig. [7] 
by P2-element 


N 

□ 

order 

0 

order 

■ 

order 

8 

7.45E-02 


4.74E-02 


1.88E-01 


16 

1.81E-02 

2.04 

6.77E-03 

2.81 

1.25E-01 

0.59 

32 

3.94E-k00 

-7.77 

1.16E-03 

2.55 

1.22E-k02 

-9.93 

64 

l.lOE-kOO 

1.84 

1.17E-04 

3.31 

8.17E-k01 

0.58 

N 


order 


order 

• 

order 

8 

6.76E-02 


4.17E-02 


1.03E-01 


16 

1.69E-02 

2.00 

8.80E-03 

2.24 

4.68E-02 

1.14 

32 

4.24E-03 

1.99 

2.04E-03 

2.11 

1.67E-02 

1.49 

64 

1.05E-03 

2.01 

4.43E-04 

2.20 

5.84E-03 

1.52 


are consistent with the theoretical ones of Scheme GSLG, Ei ^2 = 0(At + /i^ + 
and Ejji = 0{At + h'^ + h^). 

Figure [7] shows the log-log graphs of Ei ^2 and E^i versus h with v = 10“®. The 
left graph shows the results of Pi-element and Table [6] shows the values of them. 
The convergence orders of 75^2 with At = 0{h) are almost 1 in both schemes (C, 
). The orders of 75^2 with At = 0{h'^) are almost 2 for small h in both schemes 
(E, ). The orders of Ej^i are almost 1 in both schemes (■, •). The right graph 

of Fig. [7] shows the results of P2-element and Table [7] shows the values of them. 
The errors Ej ^2 with At = 0{h^) are too large at TV = 32 and 64 to be plotted in 
the graph in Scheme LG' (C) while the convergence order is almost 2 in Scheme 
GSLG ( ). The order E ]^2 with At = 0{h^) is almost 3 for small h in Scheme LG' 
(E) and almost 2 in Scheme GSLG ( ). The errors E^i are too large at A = 32 
and 64 to be plotted in the graph in Scheme LG' (■) while we can observe the 
convergence but the order is less than 2 in Scheme GSLG (•). In order to obtain 
the theoretical convergence order it seems that finer mesh will be necessary. 

7. Conclusions 

We have presented a genuinely stable Lagrange-Galerkin scheme for convection- 
diffusion problems. In the scheme locally linearized velocities are used and the 
integration is executed exactly without numerical quadrature. For the P^-element 
we have shown error estimates of 0{At + h? + in £°°(L^)-norm and of 0{At + 

-I- h^) in .^°°(iJ^)-norm. We have also obtained error estimate, c{At -\- h? + h^) 
in where the coefficient c is dependent on the exact solution (jj but 

independent of the diffusion constant v. Numerical results have reflected these 
estimates. The extension to the Navier-Stokes equations will be discussed in a 
forthcoming paper. 
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